First-principles study of improper ferroelectricity in TbMn03 
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We have carried out a first-principles theoretical study of the magnetically induced polarization 
in orthorhombic TbMnO:;, a prototypical material in which a cycloidal spin structure generates an 
electric polarization via the spin-orbit interaction. We compute both the electronic and the lattice- 
mediated contributions to the polarization and find that the latter is strongly dominant. We analyze 
the spin-orbit induced forces and lattice displacements from both atomic and mode-decomposition 
viewpoints, and show that a simple model based on nearest Mn-Mn neighbor Dzyaloshinskii-Moriya 
interactions is not able to account fully for the results. The direction and magnitude of our computed 
polarization are in good agreement with experiment. 

PACS numbers: 75.80.+q,77.80.-e 



Multifcrroic and magnetoelectric materials are cur- 
rently under intensive study, both from the point of view 
of fundamental materials physics and because of their po- 
tential application in novel technological devices such as 
memories, sensors, and transducers [lj. An interesting 
class of systems are "improper ferroelectrics" in which 
the magnetic order induces an electric polarization in 
a material that is otherwise structurally centrosymmet- 
ric. This can happen in collinear spin systems if the 
spin pattern breaks inversion symmetry; examples in- 
clude i?Mn 2 5 (R=Tb, Y, etc.) the E phase of or- 
thorhombic i?Mn0 3 (R=Ro, Tb) [J,li|/rbMn03 in mag- 
netic field [6(, and Ca 3 Co2-a ; Mn :c 06 [21] ■ In such cases, 
symmetric exchange interactions can explain the devel- 
opment of polar order, i.e., via exchange striction effects 
[a, Q . Here we shall be concerned with another class of 
systems in which a non-collinear cycloidal spin structure 
induces an electric polarization via the spin-orbit (SO) 
interaction [§]; examples include Ni 3 V 2 8 [13], CuFe0 2 

n|, and orthorhombic i?Mn0 3 (R=Tb, Dy, Ho, etc.) 

12 . |l3L Il4j . While such spiral spin structures and the 



resulting polarizations usually appear only at low tem- 
perature, these systems are of special interest because 
the coupling between magnetic and electric degrees of 
freedom is so profound. 

Theoretically, symmetry analysis and phenomenologi- 
cal models have clarified the circumstances under which 
a spiral spin structure can give rise to an electric polar- 
ization [la . Ila . Il7| , but are not well suited to identifying 
the dominant microscopic mechanism responsible for the 
polar order. Microscopic models have been introduced 
and used to discuss possible mechanisms, which can be 
divided broadly into two categories. The first are purely 
electronic mechanisms, in which the SO interaction (SOI) 
modifies the hybridization of electronic orbitals in such a 
way as to shift the center of charge 18, 111,0, El]. The 
second are "lattice mechanisms" involving magnetically- 
induced ionic displacements, usually discussed in terms 
of Dzyaloshinskii-Moriya (DM) interactions 17, 22, 23j |. 
However, it is generally difficult to estimate the magni- 



tudes (and even signs) of these two kinds of contributions 
from the models. Because experiments have not been 
sensitive enough to resolve the tiny SO-induced atomic 
displacements, a central unanswered question is whether 
electronic or lattice mechanisms are dominant. It is also 
unclear whether the pattern of displacements should fol- 
low the predictions of a simple DM model, or whether 
more complicated interactions enter the picture. 

First-principles calculations can play an invaluable role 
in resolving such questions. Density- functional theory 
(DFT) has proven to be an extremely useful tool in un- 
derstanding the physics of ordinary ferroelectric materi- 
als, and has already been used to make important con- 
tributions to the understanding of improper mag netic 
ferroelectrics of the collinear-spin type [2J, [25|, [26( and 
in the spiral magnetic materials LiCu202 and LiCuV04 
27]. However, to our knowledge such methods have not 



previously been applied to study the cycloidal magnetic 
states in the orthorhombic perovskite system. 

In this Letter, we choose orthorhombic TbMn0 3 as 
a paradigmatic system of the cycloidal type, and carry 
out a detailed first-principles study of the mechanisms by 
which the cycloidal magnetic order gives rise to the polar- 
ization. We demonstrate that DFT calculations predict 
the appearance of electric polarization in the cycloidal 
state via two mechanisms, a purely electronic one and a 
lattice-mediated one, and that the latter is strongly dom- 
inant. Nevertheless, we find that the DM theory alone 
cannot explain the computed pattern of atomic displace- 
ments, suggesting that other interactions not included in 
simple DM models play an important role. 

Ferroelectricity coupled to spiral magnetic order was 
discovered in the orthorhombic perovskite TbMn0 3 
(space group Pbnm) [12J with a 20-atom cell contain- 
ing four formula units. The Mn 3+ ions form a collinear 
sinusoidal spin wave at temperatures between ~27K and 
~41 K, but below ~27K a cycloidal spin wave forms with 
incommensurate wavevector k s ~ 0.28 along b, and a po- 
larization simultaneously appears along c [3, [l2| • 

In this work, electronic-structure calculations arc car- 
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FIG. 1: (Color online) Sketch of a x 36 x c/2 orthorhombic 
cell of TbMnOs (Pbnm space group) showing Mn 3+ sites with 
a cycloidal magnetic structure. 



TABLE I: Experimental (Ref. I35T) and theoretical (this work) 
structural parameters for orthorhombic TbMnOs. 







Experiment 


Theory 


Lattice vectors 


a (A) 


5.293 


5.195 




6(A) 


5.838 


5.758 




c(A) 


7.403 


7.308 


Tb 4c(x y 1/4) 


X 


0.983 


0.979 




V 


0.082 


0.084 


Mn 4b(l/2 0) 








01 4c(x y 1/4) 


X 


0.104 


0.107 




V 


0.467 


0.469 


02 8d(a; y z) 


X 


0.704 


0.699 




y 


0.326 


0.320 




z 


0.051 


0.052 



ried out using a plane-wave pseudopotential approach to 
DFT as implemented within the VASP code [28| using 
PAW potentials [H, H^. (The Tb potential does not 
have / electrons in the valence.) We use the local- 



density approximation (Ceperley-Alder[3l| with Vosko- 
Wilk-Nusair correlation |32|) with on-site Coulomb inter- 
actions (LDA+U). We let U=l eV to make the band gap 
close to the experimental value of ~ 0.5 eV |33[. Since 
the first-principles calculations require periodic bound- 
ary conditions, we study a system in which the cycloidal 
spin structure is commensurate with k s = 1/3, making 
a supercell that is tripled in ^-direction (60 atoms per 
cell) as shown in Fig. [I] This is close enough to the ex- 
perimental k s of 0.28 that our results can be expected to 
be qualitatively and semiquantitatively meaningful. A 
3x1x2 fc-point sampling is used. The plane- wave energy 
cutoff was 500 eV, and the electric polarization was cal- 
culated using the Berry-phase approach |34( . 

A structural relaxation was performed for the 60-atom 
supercell without SO. The resulting structural coordi- 
nates, reported in Table HI are in quite satisfactory agree- 
ment with experiment. Our computed magnetic mo- 
ments have magnitudes 3.4 /ib and rotate in the b — c 
plane as one scans along 6, and this result hardly changes 
when SO is turned on. Experimentally, the magnetic mo- 
ments on Mn 3+ sites are found to form an elliptical spiral 
with amplitudes m& = 3.9 /^b and m c = 2.8 hb The 
origin of this ellipticity, which is apparently not captured 
by our calculation, deserves further investigation. 

Next, we confirmed that the Berry-phase polarization 
vanishes in the absence of SO. When the SO is turned 
on, one expects the broken inversion symmetry in the 
spin sector to be communicated to the spatial (charge) 
degrees of freedom, and indeed our calculations confirm 
this. Focusing first on the purely electronic mechanism, 
we froze the lattice vectors and ionic positions while turn- 
ing on the SOI and computed the polarization to be 
P = 32^C/m 2 along the c direction[36j (and zero along 
directions a and 6) . This is much smaller than the exper- 
imental value of ~ 600/iC/m 2 [12|, suggesting that the 



lattice mechanism must be dominant. We thus turned 
our attention next to the role of ionic displacements. 

As a first step, we computed the Hellmann-Feynman 
forces appearing on the ions as the SO is turned on 
(with the structural coordinates still clamped). Using 
the 20-atom cell of the Pbnm space group to label the 
forces, we found that: (i) 64.8% of the forces belong to 
the symmetry-preserving Ai(r) irreducible representa- 
tion (irrep); (ii) 28.7% of the forces belong to the B3 U (r) 
irrep, i.e., zone-center infrared-active (IR-active) modes 
with dynamical dipoles along c; and (iii) 6.5% correspond 
to forces at k y = ±2ir/3b, i.e., those that would be re- 
sponsible for lowering the translational symmetry of the 
20-atom cell. The remaining contributions are found to 
vanish to numerical precision. Since the only forces that 
generate an electric polarization in linear order are those 
of type (ii), we henceforth focus our attention on the 
zone-center IR-active B3 U modes. 

Eight Wyckoff coordinates contribute to the B3 U irrep 
in the 20-atom Pbnm structure: three associated with Mn 
atoms in Wyckoff position 4b, one each for Tb and Ol 
atoms in Wyckoff position 4c, and three for the 02 atoms 
in Wyckoff position 8d jsjj]. To find the polarization in- 
duced by the SOI, we need to deduce the displacements 
resulting from these forces. To this end we calculated 



the 8x8 force-constant matrix 



-dfi/duj connect- 



ing the B3 U forces and displacements at linear order using 
finite-difference methods (i.e., we shift each Wyckoff co- 
ordinate by ^10 -3 A and calculate the resulting forces). 
For computational convenience, the determination of $jj 
was done in the absence of SO. (We checked a few el- 
ements and found hardly any difference in the presence 
of SO 38].) Furthermore, we identify the linear com- 
bination tj of the eight coordinates corresponding to a 
uniform shift of the crystal along the c direction (acous- 
tic mode) and make sure to project this combination out 
of the computed forces and displacements. Thus, the 
displacements are predicted from Ui = — Y2j=i ^Tj'fi' 
where $ _1 is a pseudo-inverse (39| having the properties 
that = u for any u orthogonal to t. 
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TABLE II: Forces, displacements, effective charges and con- 
tributions to ionic polarization from IR-active E>3 U modes. 



Wyckoff 


F 


Au 


Z* 


AP 


coordinate 


(meV/A) 


(ml) 


(e) 


( M C/m 2 ) 


1 (Tb 4c, z) 


0.434 


-0.171 


7.468 


-94 


2 (01 4c, z) 


2.264 


-0.162 


-6.815 


81 


3 (Mn 4b, x) 


-7.039 


-0.315 


0.567 


-13 


4 (Mn 4b, z) 


-8.933 


-0.453 


7.460 


-248 


5 (Mn 4b, y) 


-2.940 


-0.012 


0.553 


-0.5 


6 (02 8d, x) 


5.058 


0.540 


0.077 


3 


7 (02 8d, y) 


3.573 


0.934 


0.227 


16 


8 (02 8d, z) 


4.409 


0.556 


-5.738 


-234 



TABLE III: Forces, displacements, effective charges and con- 
tributions to ionic polarization from modes corresponding to 
eigenvectors of the force-constant matrix. 



Eigenvalue 


F 


Au 


Z* 


AP 


(eV/A 2 ) 


(meV/A) 


(mA) 


(e) 


(pC/m 2 ) 


34.86 


6.349 


0.182 


-8.111 


-108 


31.37 


-8.259 


-0.263 


4.565 


-88 


13.32 


2.783 


0.209 


8.132 


125 


11.87 


7.614 


0.641 


-2.127 


-100 


10.35 


1.684 


0.163 


-3.406 


-41 


6.99 


-2.784 


-0.398 


4.316 


-126 


3.85 


4.032 


1.047 


-1.991 


-153 



The computed forces and predicted atomic displace- 
ments in the IR-active B3 U sector for TbMn0 3 are pre- 
sented in the first three columns of Table [TTJ We calcu- 
lated the B 3u (r) forces again in the presence of the pre- 
dicted displacements and found that 99% of these forces 
were eliminated, thus justifying our approach (40| . The 
distorted crystal structure is now non-centrosymmetric 
(space group Pna2\), and the calculated Berry-phase po- 
larization (including both electronic and ionic contribu- 
tions) is P — —467 /iC/m 2 . This is more than an order of 
magnitude larger than the value of 32/^C/m 2 computed 
from purely electronic mechanisms, and comparable in 
magnitude with the experimental value of ~ —600 /iC/m 2 
Note that the sign of the ionic contribution is dif- 
ferent from the purely electronic contribution, making 
the sign of the total polarization consistent with exper- 
iment. We thus arrive at the first major conclusion of 
our paper, namely, that the lattice contributions domi- 
nate strongly and have approximately the right magni- 
tude to explain the experimentally observed polarization. 
Purely electronic mechanisms, e.g. as described by the 
Katsura-Nagaosa-Balatsky (KNB) model[I3l or related 



work [H, Efl E3, are not sufficient to describe the ferro- 
electricity in TbMn0 3 . 

We focus henceforth on the lattice-mediated contribu- 
tion to the polarization. To understand which modes 
contribute most strongly and to test whether simple mod- 
els (based, e.g., on DM interactions) can predict the 
forces and displacements, we performed a further anal- 
ysis as follows. For each B 3u Wyckoff coordinate, we 
computed the corresponding effective charge Z* (deriva- 
tive of polarization with respect to the Wyckoff coordi- 
nate) by augmenting the same finite-difference calcula- 
tions used earlier with Berry-phase polarization calcu- 
lations. We then multiplied effective charges and dis- 
placements to get the contribution AP made by each 
Wyckoff coordinate. The results are given in the last two 
columns of Table [TXJ The total ionic polarization (sum 
of the last column) is — 489 ^C/m 2 , in good agreement 
with the change of polarization in going from the cen- 
trosymmetric to the distorted structure as reported ear- 
lier (-467-32=-499^C/m 2 ). Not surprisingly, the Z* 



values are much larger for the four Wyckoff coordinates 
involving displacements along z, and the corresponding 
AP values give by far the largest contribution. 

In a simple DM model [f|, one considers interactions 
involving nearest-neighbor triplets of ions along Mn-O- 
Mn bonds. Such a picture leads to the expectation that 
a force 7 e nn r x S„ x S„; should appear on the central O 
ion (where 7 is a coupling, S„ and S ra ' are neighboring 
Mn spins, and e nn i is the unit vector connecting them), 
with an equal and opposite force shared between the Mn 
ions. The cycloidal spin structure of TbMn0 3 is such 
that 01 ions should be unaffected, and 02 ions should 
all feel the same force along c, with Mn ions feeling an 
opposite force. Thus, within such a model one would 
expect only the 4th and 8th rows of Table [III to give 
significant contributions. Inspecting the Table, we see 
that this is rather far from being the case; one does find 
substantial oppositely-directed forces and displacements 
for the (Mn 4b, z) and (02 8d, z) coordinates, but the 
forces and displacements are comparable for some other 
coordinates. The contributions to AP are more strongly 
dominated by these two coordinates, but only because 
the Z* values are so much larger for displacements along 
z. Interestingly, the Tb displacements are responsible 
for ^20% of the lattice polarization, despite being unin- 
volved in the usual picture of Mn-O-Mn couplings. Thus 
we conclude that the simple DM model, and indeed any 
model that considers only nearest-neighbor Mn-Mn spin 
interactions, is unable to provide a detailed description 
of the ferroelectricity in TbMn0 3 . 

As an alternative way of analyzing the SO-induced B 3u 
forces and displacements, we transform into the basis of 
eigenvectors of the 8x8 force-constant matrix <E>. The 
results are presented in Table IIIII (Components on the 
translational zero mode, not shown in the Table, were nu- 
merically negligible.) We find that the forces are rather 
widely distributed among these "modes." Not surpris- 
ingly, the displacements are largest for the softest mode, 
but other modes also show substantial displacements. 
Moreover, the mode dynamical charges are much larger 
for some of the harder modes, with the result that the 
contributions to the polarization are of comparable mag- 
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nitude for almost all of the B3 U modes. Thus, it is clear 
that the physics of the electric polarization in TbMnO,3 
is not dominated by a single softest mode. 

To investigate the role of spin-orbit effects located on 
the various atomic sites, we performed additional calcu- 
lations of the forces (and, via the force-constant matrix, 
the resulting lattice polarizations). We first turned the 
SOI off on all sites other than Tb, then Mn, then O, 
finding contributions to the polarization of — 11/LtC/m 2 , 
— 447 ^C/m 2 and — 8/iC/m 2 , respectively. (The sum is 
in fairly good agreement with the total ionic polarization 
of — 489/iC/m 2 ; this and other tests confirm that the re- 
sults are roughly linear in SO coupling.) This confirms 
that the SOI on the Mn sites is responsible for almost all 
of the lattice-mediated contribution to the polarization, 
even though some of these forces appear on other sites, 
even including Tb (see Table HT|) . 

Up to now our analysis has been done using the 60- 
atom supercell, whose cycloidal wavevector (fc s =l/3) ap- 
proximates the experimental one (fc s =0.28). To investi- 
gate how the results depend on k s , we have also carried 
calculations of the SO-induced forces for a 40-atom super- 
cell (fc s =l/2). In order to isolate the effects of the spin- 
wave period, we calculated the forces at the same stru- 
cural coordinates as for the 60-atom structure, and again 
projected into the B 3u (r) irrep. We find that the pattern 
of B3 U forces matches almost exactly that of the 60-atom 
cell (with the angle between the two eight-dimensional 
vectors being ~3°), although the magnitude is larger for 
the 40-atom cell. From the simple DM model one ex- 
pects P oc e n>n i x (S„ x S„/) [9J where the spins are 
neighbors on the zig-zag chains propagating along b (see 
Fig. [1]) . The angle between the spins is 90° and 60° for 
the 40-atom and 60-atom supercells respectively, giving 
an expected ratio of 2 : y/3=1.15. In the long-wavelength 
limit one expects P oc |fe s |, for a ratio of 1.5. Strangely, 
our computed ratio for the 40- and 60-atom supercells is 
about 1.48, closer to the continuum expression in spite 
of the rather short wavelength of the spiral. This may be 
another indication that further-neighbor Mn-Mn interac- 
tions, or other effects not captured by a nearest-neighbor 
DM model, play an important role. From the forces we 
predict a lattice polarization of — 754 ^(C/m 2 for fe s =l/2; 
extrapolating through our fc s =l/3 result, we estimate a 
lattice polarization of —410 /iC/m 2 for fc s =0.28. 

In summary, we have used first-principles methods to 
compute the electronic and lattice contributions to the 
spin-orbit induced electric polarization in the cycloidal- 
spin compound TbMn03, and find the lattice contribu- 
tion to be strongly dominant. We provide a detailed anal- 
ysis of the lattice contributions coming from individual 
sites and individual modes and roughly characterize the 
dependence of the lattice polarization on spin-wave pe- 
riod. Our result for the polarization agrees in sign and 
compares fairly well in magnitude with the experimen- 
tally measured value of ~ — 600/iC/m 2 14 1. 



The authors wish to thank Morrel H. Cohen and Craig 
J. Fennie for fruitful discussions. This work was sup- 
ported by NSF Grant No. DMR-0549198. 

Note: While finalizing this manuscript, we learned of 
closely related work by H.J. Xiang et al. [4l[ 
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